Boundary multifractality in critical ID systems with long-range hopping 



o 
o 



o 



A. Mildenberger 1 , A. R. Subramaniam 2 , R. Narayanan 3 ' 5 , F. Evers 3 ' 4 , I. A. Gruzberg 2 , and A.D. Mirlin 3 ' 4 '* 
1 Fakultat fur Physik, Universitdt Karlsruhe, 76128 Karlsruhe, Germany 
2 James Franck Institute, University of Chicago, 5640 S. Ellis Ave., Chicago, IL 60637 
3 Institut fur Nanotechnologie, Forschungszentrum Karlsruhe, 76021 Karlsruhe, Germany 
4 Inst, fur Theorie der Kondensierten Materie, Universitdt Karlsruhe, 76128 Karlsruhe, Germany 
5 Department of Physics, Indian Institute of Technology Madras, Chennai 600036, India 

(Dated: February 3, 2008) 

Boundary multifractality of electronic wave functions is studied analytically and numerically for 
the power-law random banded matrix (PRBM) model, describing a critical one-dimensional system 
with long-range hopping. The peculiarity of the Anderson localization transition in this model is the 
existence of a line of fixed points describing the critical system in the bulk. We demonstrate that 
the boundary critical theory of the PRBM model is not uniquely determined by the bulk properties. 
Instead, the boundary criticality is controlled by an additional parameter characterizing the hopping 
amplitudes of particles reflected by the boundary. 
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I. INTRODUCTION 

Although almost half a century has passed after An- 
derson's seminal paper—, the properties of disordered sys- 
tems at Anderson localization transitions remain a sub- 
ject of active current research. This interest is addition- 
ally motivated by the understanding that, besides the 
"conventional" Anderson transition in d > 2 dimensions, 
disordered fermions in two dimensions (2D) possess a rich 
variety of critical points governing quantum phase tran- 
sitions in these systems. 

One of the striking peculiarities of the Anderson 
transitions is the multifractality of the electronic wave 
functions 2 -! 3 -, see Refs. |4||5| for recent reviews. Specifically, 
the scaling of moments of the wave functions with the 
system size L is characterized by a continuum of inde- 
pendent critical exponents r q , 



^(l^(i-)l 29 ) 



Tq = d(q - 1) + A„ 



(1) 



where (. . .) denotes the disorder average. In the field- 
theoretical (c-model) language^ r q are scaling dimen- 
sions of higher-order operators describing the moments of 
the local density of states. Note that one often introduces 
fractal dimensions D q via r q = D q (q — 1). In a metal 
D q = d, while at a critical point D q is a nontrivial func- 
tion of q, implying the multifractality of wave functions. 
Nonvanishing anomalous dimensions A q = (q—\)(D q — d) 
distinguish a critical point from a metallic phase and de- 
termine the scaling of wave function correlations. Among 
them, A2 < plays the most prominent role, governing 
the spatial correlations of the intensity \ip\ 2 , 



L 2d (|^(r)V 2 (r')|)^(|r-r'|/i) Z 



(2) 



This equation, which in technical terms results from 
an operator product expansion of the field theory^, can 
be obtained from ([I} by using the fact that the wave 
function amplitudes become essentially uncorrelated at 
|r — r'| ~ L. Scaling behavior of higher order spatial cor- 
relations, (\yb 2qi (ri)ip 2q2 (r 2 ) . . . yb 2qn (r n )\), can be found 



in a similar way. Above, the points rj were assumed to 
lie in the bulk of a critical system. In this case we denote 
the multifractal exponents by T q , A q , etc. The multifrac- 
tality of wave functions has been studied analytically and 
numerically for a variety of systems: Anderson transition 
in d = 2 + e, 3, and 4 dimensions 2 -^, as well as weak 
multifractality 9 , Dirac fermions in random gauge fields^, 
symplectic-class Anderson transition^, integer quantum 
Hall 12 and spin quantum Hall 13 transitions in two dimen- 
sions, and power-law random banded matrices". 

Recently, the concept of the wave function multifrac- 
tality was extended^ to the surface of a system at the 
critical point of an Anderson transition. It was shown 
that the fluctuations of critical wave functions at the 
surface are characterized by a new set of exponent 
(or, equivalently, anomalous exponents A®), which are in 
general independent from their bulk counterparts. This 
boundary critical behavior was explicitly studied, ana- 
lytically as well as numerically, for the 2D spin quantum 
Hall transition 1 - 5 ^ and a 2D weakly localized metali 5 - 
and, most recently, for the Anderson transition in a 2D 
system with spin-orbit couplingi 7 -. 

In the present paper, we analyze the boundary criti- 
cality in the framework of the power-law random banded 
matrix (PRBM) model 5 -^^. The model is defined^ as 
the ensemble of random Hermitian matrices H (real for 
(3 = 1 or complex for (3 = 2). The matrix elements Hij 
are independently distributed Gaussian variables with 



zero mean 



(Hi. 



and the variance 



a 2 (\i-j\) 



where a(r) is given by 

a 2 (r) = 



1 + (r/b) 



2a 



(3) 



(4) 



At a = 1 the model undergoes an Anderson transition 
from the localized (a > 1) to the delocalized (a < 1) 
phase. We concentrate below on the critical value a = 1 , 
when a(r) falls down as a(r) a 1/r at r > 6. 
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In a straightforward interpretation, the PRBM model 
describes a ID sample with random long-range hopping, 
the hopping amplitude decaying as l/r Q with the dis- 
tance. Also, such an ensemble arises as an effective de- 
scription in a number of physical contexts (see Ref. [lj 
for relevant references). At a = 1 the PRBM model 
is critical for arbitrary values of b and shows all the 
key features of an Anderson critical point, including 
multifractality of eigenfunctions and nontrivial spectral 
compressibility 5,18 . The existence of the parameter b 
which labels critical points is a distinct feature of the 
PRBM model: Eq. © defines a whole family of critical 
theories parametrized by b. The limit b 1 represents 
a regime of weak multifractality, analogous to the con- 
ventional Anderson transition in d = 2 + e with e -C 1. 
This limit allows for a systematic analytical treatment via 
the mapping onto a supermatrix er-model and the weak- 
coupling expansio n 5 ' 14 ! 18 . The opposite limit b <C 1 is 
characterized by very strongly fluctuating eigenfunctions, 
similarly to the Anderson transition in d ^S> 1, where the 
transition takes place in the strong disorder (strong cou- 
pling in the field-theoretical language) regime. It is also 
accessible to an analytical treatment using a real-space 
renormalization-group (RG) method 14 introduced earlier 
for related models in Ref. [191 

In addition to the feasibility of the systematic ana- 
lytical treatment of both the weak-coupling and strong- 
coupling regimes, the PRBM model is very well suited 
for direct numerical simulations in a broad range of cou- 
plings. For these reasons, it has attracted a consider- 
able interest in the last few years as a model for the 
investigation of various properties of the Anderson criti- 
cal poin t 14 ' 20 . We thus employ the PRBM model for the 
analysis of the boundary multifractality in this work. The 
existence of a line of fixed points describing the critical 
system in the bulk makes this problem particularly inter- 
esting. We will demonstrate that the boundary critical 
theory of the PRBM model is not uniquely determined by 
the bulk properties. Instead, the boundary criticality is 
controlled by an additional parameter characterizing the 
hopping amplitudes of particles reflected by the bound- 
ary. 

The structure of the paper is as follows. In Sec. HT1 we 
formulate the model. Section ITTT1 is devoted to the ana- 
lytical study of the boundary multifractal spectrum, with 
the two limits b 3> 1 and b <C 1 considered in Sees. IIII Al 
and IIII Bl respectively. The results of numerical simula- 
tions are presented in Scc. lIVI Section fVl summarizes our 
findings. 



II. MODEL 

We consider now the critical PRBM model with a 
boundary at i = 0, which means that the matrix ele- 
ment Hij is zero whenever one of the indices is negative. 
The important point is that, for a given value of the bulk 
parameter b, the implementation of the boundary is not 



unique, and that this degree of freedom will affect the 
boundary criticality. Specifically, we should specify what 
happens with a particle which "attempts to hop" from a 
site i > to a site j < 0, which is not allowed due to the 
boundary. One possibility is that such hops are simply 
discarded, so that the matrix element variance is simply 
given by {\H tJ | 2 ) = [1 + (i - j) 2 /b 2 }- 1 for i, j > 0. More 
generally, the particle may be reflected by the boundary 
with certain probability p and "land" on the site — j > 0. 
This leads us to the following formulation of the model: 



Jij — 



1 



-3\ 2 /b 2 



+ 



P 



i + V + jV /b 2 



(5) 
(6) 



While the above probability interpretation restricts p to 
the interval [0, 1], the model is defined for all p in the 
range — 1 < p < oo. The newly introduced parameter p is 
immaterial in the bulk, where i, j 3> \i— j\ and the second 
term in Eq. © can be neglected. Therefore, the bulk 
exponents depend on b only (and not on p), and their 
analysis performed in Ref. remains applicable without 
changes. On the other hand, as we show below by both 
analytical and numerical means, the surface exponents 
Tg are a function of two parameters, b and p. 

Equation ([6]) describes a semi-infinite system with one 
boundary at i = 0. For a finite system of a length L (im- 
plying that the relevant coordinates are restricted to < 
i,j < L) another boundary term, p'/[l + (i+j — 2L) 2 /b 2 ], 
is to be included on the right-hand side of Eq. ([6]). In 
general, the parameter p' of this term may be different 
from p. This term, however, will not affect the boundary 
criticality at the i — boundary, so we discard it below. 



III. BOUNDARY MUTLIFRACTALITY: 
ANALYTICAL METHODS 

A. 6>1 

The regime of weak criticality, b 3> 1, can be stud- 
ied via a mapping onto the supermatrix g-model- i 14 ' 18 , 
in analogy with the conventional random banded matrix 
model^i. The cr-model action has the form 



S[Q] = ^Str 



{■nv) 2 2J JjjQiQj - ittuoj QjA 

i,j=0 i=0 



(7) 

where Q r is a 4 x 4 ((3 — 2) or 8 x 8 ((3 = 1) supermatrix 
field constrained by Q 2 — 1, A = diag(l,— 1), and Str 
denotes the supertracei^ 2 Furthermore, J,j are given by 
Eq. ©, u> is the frequency, and v is the density of states 
given by the Wigner semicircle law 

V{E)= 2^b {47rb ~ E2)1/2 ' 1^1 < ( 8 ) 

For definiteness, we will restrict ourselves to the band 
center, E = 0, below. 
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To calculate the multifractal spectrum to the leading 
order in 1/6 <C 1, we will need the quadratic form of the 
action expressed in terms of independent coordinates. 
Parametrizing the field Q (constrained to Q 2 = 1) in the 
usual way, 



Qi 



A 1 + Wi 



2 



(9) 



we obtain the action to the second order in the W fields, 



S[W] 



where 



-Str Wi ^(J^Sij-Jij) 



r(i) 



k=0 



(10) 

(11) 



The equation of motion for this action reads (after the 
Fourier transformation from the frequency into the time 
domain) 



dWjjt) 
dt 



\S- ■ 7° 



W j (t) = 0. 



(12) 



This equation is the analog of the diffusion equation for 
a metallic system. 

The cr-model action allows us to calculate the moments 
(\ip 2 \ q ) at a given point r. On the perturbative level, the 
result reads^ii^i2i 



(ivv 2 n = m 



2|<?\ 



'RMT 



1 



1 + pQ(q ~ l)H r 



(13) 



Here the factor (|VvI 9 )rmt ^ s ^ ne random-matrix-theory 
result equal to (2q — l)UL~ q for (3 — 1 and q\L~ q for f3 = 
2. The second term in the square brackets in Eq. (fT5)l . 
which constitutes the leading perturbative correction, is 
governed by the return probability Tl rr to the point r, i.e., 
the diagonal matrix element of the generalized diffusion 
propagator Tl rr i . The latter is obtained by the inversion 
of the kinetic operator of Eqs. (flU)) . (TH 



oo 

■£[ 

3=0 



6--J {t} 



Jij 



n 



(14) 



(The "diffusion" operator has a zero mode related to the 
particle conservation. The term L _1 on the right-hand 
side of Eq. JT3]) ensures that the inversion is taken on 
the subspace of nonzero modes.) In the bulk case, the 
inversion is easily performed via the Fourier transform, 



n,. 



with t 



-(nv) 2 b 2 , 



t 



\k\ < b-\ 



(15) 



A/b at the band cen- 



ter. The l/\k\ behavior of the propagator should be 



contrasted to its 1/fc 2 scaling for a conventional metallic 
(diffusive) system. This implies that the kinetics gov- 
erned by Eq. (|12p is superdiffusive, also known as Levy 
flights^. Substitution of JTSJ) in Eq. HHJ yields a loga- 
rithmic correction to the moments of the wave function 
amplitude, 



2 I?\ 



(HIV) 



RMT 



1 , g(g-l) , L 
1 + ^3T ln I 



(16) 



Equation (fT6|l is valid as long as the relative correction is 
small. The logarithmic divergence of the return probabil- 
ity in the limit L — > oo, which is a signature of criticality, 
makes the perturbative calculation insufficient for large 
enough L. The problem can be solved then by using the 
renormalization group (RG)r ' 14 i 18 which leads to the ex- 
ponentiation of the perturbative correction in Eq. (|16j) . 
This results in Eq. (fl]) with the bulk multifractal expo- 
nents 



(17) 



(l-l) 1 



2tt/36 



The first term (unity) in the second factor in Eq. (fT7|) 
corresponds to the normal (metallic) scaling, the second 
one determines the anomalous exponents 



b _ 9(1-9) 



A 9 = 



2tt/% 



(18) 



At the boundary, the behavior is qualitatively the same: 
the return probability Tl rr increases logarithmically with 
the system size L, in view of criticality. However, as we 
show below, the corresponding prefactor [and thus the 
prefactor in front of the second term in square brackets 
in Eq. (fT6|) ] is different. After the application of the RG 
this prefactor emerges in the anomalous exponent, 



g(l-g) ff 
27T/36 p 



A h q R p . 



(19) 



In the presence of a boundary the system is not transla- 
tionally invariant anymore, which poses an obstacle for 
an analytical calculation of the return probability H rr . 
While for Levy-flight models with absorbing boundary 
(that is obtained from our equation (|12jl with p = by a 
replacement of with its bulk value J ) an analytical 
progress can be achieved via the Wiener-Hopf method 24 , 
it is not applicable in the present case, since the kernel 
of Eq. (Tj"2"|) is not a function of i — j only. We thus pro- 
ceeded by solving the classical evolution equation JT2]) 
numerically with the initial condition Wi(0) = 8; Lr . The 
value W r {t) of the solution at the point r (i.e. the prob- 
ability to find the particle at the initial point) decays 
with the time as 1/i, so that the integral J dtW r (t) yields 
the logarithmically divergent return probability discussed 
above. Extracting the corresponding prefactor, we find 
the anomalous exponent, 



(20) 
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FIG. 1: The ratio R v = A% (b, p) /A° (6) of the surface and bulk 
anomalous exponents for large b, as a function of the reflection 
parameter p. Diamonds represent the results of the cr-model 
analysis with a numerical solution of the corresponding clas- 
sical evolution equation, as described in Sec. IIII Al Circles 
represent a direct computer simulation of the PRBM model, 
Eqs. (J5]l, ©, see Sec. [IV] with 6 = 8. The ratio R p has been 
evaluated for the range < q < 1, where the numerical ac- 
curacy of the anomalous exponents is the best. Within this 
interval we find that R p is g-independent (within numerical 
errors) in agreement with Eq. 1)190 , 



Note that the limit of the large system size L — > oo should 
be taken in Eq. (j2"0|) before t — > oo, so that the particle 
does not reach the boundary for r in the bulk (or, for r 
at the boundary, does not reach the opposite boundary) . 

We have checked that the numerical implementation of 
Eqs. (|20|) . (fT2"l) reproduces the analytical result (fT5|) in the 
bulk. We then proceeded with numerical evaluation of 
the surface multifractal exponents A^. For this purpose, 
we have discretized the time variable in Eq. (|12p with a 
step At = 1/2. With the parameters itv = 1/ir, b = 10, 
L = 10000, and t = 500, the product tW r (t) yields its 
required asymptotic value with the accuracy of the order 
of 2%. The results for the corresponding prefactor R p , 
as defined in Eq. (fl9|) , are shown in Fig. Q] for several 
values of p between and 3. It is seen that the boundary 
exponents not only differ from their bulk counterparts 
but also depend on p. 

For the particular case of the reflection probability p — 
1 we can solve the evolution equation (fT2"|) and find A^ 
analytically. Indeed, the corresponding equation can be 
obtained from its bulk counterpart (defined on the whole 
axis, — oo < i < oo) by "folding the system" on the 
semiaxis i > according to Wi{t) + W-i(t) — > Wi(t), 
cf. Ref. 25. This clearly leads to a doubling of the return 
probability, so that 



i?i =2, 



(21) 



in full agreement with the numerical solution of the evo- 
lution equation. 



In the regime of small b the eigenstates are very sparse. 
In this situation, the problem can be studied by a real- 
space RG method that was developed in Ref. [l9| for re- 
lated models and in Ref. [13 for the PRBM model. Within 
this approach, one starts from the diagonal part of the 
Hamiltonian and then consecutively includes into consid- 
eration nondiagonal matrix elements with increasing 
distance p = \i — j\. The central idea is that only rare 
events of resonances between pairs of remote states are 
important, and that there is an exponential hierarchy of 
scales at which any given state finds a resonance part- 
ner: lnpi ~ b" 1 , ln^2 ~ b~ 2 , .... This allows one to 
formulate RG equations for evolution of quantities of in- 
terest with the "RG time" t = In p. We refer the reader 
for technical details of the derivation to Ref. where 
the evolution equation of the distribution f(P q ) of the 
inverse participation ratios, P q — J2 r \^r\ q ^ as weu as 
of the energy level correlation function, was derived. In 
the present case, we are interested in the statistics of the 
local quantity, the wave function intensity \t/j^\ at a cer- 
tain point r. Assuming first that r is in the bulk and 
generalizing the derivation of Ref. [TH, we get the evolu- 
tion equation for the corresponding distribution function 

f(v=m), 



df(y,p) 
<9 hip 



2b 

71" Jn 



•72 



dO 



sin 2 9 cos 2 I 



dy'f(y',p) 



x [6(y ~ y' cos 2 8) + 6(y - y' sin 2 6) 

-5(y-y')-5(y)}. (22) 

Equation ([22|) is written for fi = 1; in the case of (3 = 2 
one should make a replacement b — > (n/2^/2)b. The 
physical meaning of Eq. (f2"2")l is rather transparent: its 
right-hand side is a "collision integral" describing a res- 
onant mixture of two states with the intensities y' and 
at the point r, leading to formation of superposition 
states with the intensities y' cos 2 9 and y' sin 2 9. Multi- 
plying Eq. (|22|) by y q and integrating over y, we get the 
evolution equation for the moments (y q ), 



where 
T(q) 



d(y q ) _ 
dlnp 

1 r /2 d9 



7T Jq sin 2 9 cos 2 
1 r(2g-l) 
22?-3 T(q)T(q - 1) 



2bT(q)(y q ) , 



-(1- cos 2<? 9 -sin 29 1 



(23) 



(24) 



The RG should be run until p reaches the system size L. 
Thus, the bulk multifractal exponents are equal to 



r\ = 2bT(q), 



(25) 



in agreement with Ref. \\A 
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How will the evolution equation (22)) be modified if the 
point r is located at the boundary? First, the factor 2 
on the right-hand side of (f22|) will be absent. Indeed, 
this factor originated from the probability to encounter 
a resonance. In the bulk, the resonance partner can be 
found either to the left or to the right, thus the factor 
of two. For a state at the boundary only one of these 
possibilities remains, so this factor is absent. Second, 
one should now take into account also the second term 
in the variance Jy of the matrix element H i j 1 Eq. JS). In 
view of the hierarchy of resonances described above, the 
relevant matrix elements will always connect two points, 
one of which is much closer to the boundary than the 
other (say, i<j). In this situation, the two terms in (6} 
become equivalent (up to the prefactor p in the second 
term) and can be combined, 



0.1 



10 



(l+p)b 2 
J 2 



i<3- 



(26) 



Therefore, the effect of the second term amounts to the 
rescaling b — ► (1 +p) 1 / 2 6. Combining both the effects, 
we get the boundary multifractal exponents, 



T s = (l+p) 1 ' 2 bT{q) = 



(1+p) 



1/2 



(27) 



The above real-space RG method works for q > 1/2, 
where the multifractal exponent r q is small2£. The results 
can, however be extended to the range of q < 1/2 by 
using the recently found symmetry relation between the 
multifractal exponents^ 7 -, 



A„ 



(28) 



Independently of whether q is larger or smaller than 1/2, 
the obtained relation between the surface and the bulk 
multifractal spectra can be formulated in the following 
way: 



T» q (b,p)=rZ(b-^ 6(1 + P ) 1/2 /2) 



IV. 



BOUNDARY MULTIFRACTALITY: 
NUMERICAL SIMULATIONS 



(29) 



In this section we present the results for the multi- 
fractality spectra obtained by direct numerical simula- 
tions of the PRBM model, Eqs. ©, ©, with = 1. 
The model has been implemented with two boundaries, 
each one having the same boundary parameter p. Us- 
ing standard diagonalization routines, systems with sizes 
L = 128, 256, 512, 1024, 2048 and 4096 sites have been 
studied with ensembles comprising of 10 7 (L — 128, 3T0 8 
wave functions) to 5000 [L = 4096, 5 • 10 6 wave func- 
tions) disorder configurations. The multifractal analysis 
has been performed with intensities \ip 2 \ averaged (coarse 
grained) over blocks of four neighboring sites in order to 
access negative q values, q > —2. For the analysis of the 
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FIG. 2: Upper panel: Anomalous exponent A2 = D2 — 1 as 
a function of b from numerical simulations in the bulk and at 
the boundary for the reflection parameter p = and 1. The 
inset shows data for p — 3 compared to the p = bulk values. 
Lower panel: Surface and bulk data for p = compared with 
analytical results for small and large b (using Ro = 2.78), 
Eqs. Jig}, (19), (25), (27). 



surface multifractal exponents, only the four sites closest 
to boundaries have been taken into account. 

Figure [5] illustrates nicely our main findings. We 
show there the dependence of the anomalous dimension 
A2 = D2 — 1 on b in the bulk and at the boundary, for 
three different values of the reflection parameter p. It 
is seen, first of all, that the bulk exponent does not 
depend on p, in agreement with the theory. Second, the 
boundary exponent A2 is different from the bulk one. 
Third, the boundary exponent is not determined by b 
only, but rather depends on the boundary parameter p as 
well. The lower panel of Fig. [2] demonstrates the agree- 
ment between the numerical results and the analytical 
asymptotics of small and large b. 

Having discussed the 6-dependence of the fractal expo- 
nent with fixed q (equal to 2) shown in Fig. [51 we turn to 
Fig.[3l where the whole multifractal spectra A 9 are shown 
for fixed large values of b. Specifically, the anomalous di- 
mensions A^ and A^ are presented for b — 2, 4, and 8, 
with the reflection parameter chosen to be p = 1. For 
all curves the q dependence is approximately parabolic, 
as predicted by the large-6 theory, Eqs. (T5|) and (TTJj) . 
with the prefactor inversely proportional to b. To clearly 
demonstrate this, we plot in the lower two panels the 
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FIG. 3: Upper panel: Boundary and bulk multifractal spec- 
tra, A B q and A q , at b = 2, 4, and 8 for the reflection parameter 
p = 1. In accordance with Eq. (|21|l . the surface multifractal- 
ity spectrum is enhanced by a factor close to two compared 
to the bulk. 

Middle panel: Surface spectrum divided by the analytical 
large-6 result, Eq. (|19|1 . The dashed line represents the ana- 
lytical result for 6 > 1. With increasing b, the numerical data 
nicely converges towards the analytical result. 
Lower panel: Analogous plot for the bulk spectrum, Eq. (|18[1 . 
The error estimate from the finite size extrapolation is 3%. 



exponents A q divided by the corresponding analytical 
results of the large-6 limit. While for moderately large 
6 the ratio shows some curvature, the latter disappears 
with increasing 6 and the ratio approaches unity, thus 
demonstrating the full agreement between the numerical 
simulations and the analytical predictions. It is also seen 
in Fig. [3] that the bulk multifractality spectrum for b = 4 
and the surface spectrum for 6 = 8 are almost identical, 
in agreement with Eq. ([2"Tj) . The same is true for the 
relation between the bulk spectrum for 6 = 2 and the 
surface spectrum for 6 = 4. 

We have further calculated the ratio of the large-6 sur- 
face and bulk anomalous dimensions, R p — A s q /A^, for 
several values of the reflection parameter, p = 0, 1, and 
3. As shown in Fig. [TJ the results are in good agree- 
ment with the cr-model predictions for R p obtained in 

sec.nrng 

In Fig. [3] the surface and bulk multifractal spectra are 
shown for the case of small 6. While the spectra are 
strongly nonparabolic in this limit, they clearly exhibit 
the symmetry q — > 1 — q, Eq. ([28]) . The data are in 
good agreement with the RG results of Sec. IIII Bl In 
particular, the surface spectrum for p — 3 is essentially 
identical to the bulk spectrum, as predicted by Eq. 
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FIG. 4: Main panel: Numerically determined boundary and 
bulk anomalous dimensions A q at b = 0.1 for p = 0, 1, and 3. 
As expected, the bulk anomalous dimension is independent of 
the value of p. In accordance with Eq. (|29[1 , for p — 3 surface 
and bulk dimensions have the same values. 
Inset: The p — data compared to the analytical results, 
surface [solid line, Eq. p7|l] and bulk [dashed line, Eq. (f25|) ]. 
Analytical data have been calculated for q > 0.6 and mirrored 
for q < 0.4 by using the symmetry relation A q — Ai_ q . Note 
that the analytical result (|25f) breaks down in the vicinity of 
q = 1/2, at |g— 1/2 < 1/ In ft -1 , see Ref.[H At the boundary, 
one should replace b — > 6(1 +p) 1 ^ 2 /2 in this condition. We 
see, indeed, that for b — 0.1, p = the analytical formula 
works up to \q - 1/2 ~ 0.4 in the bulk and \q - 1/2| ~ 0.3 at 
the boundary. 



In the inset, the surface and bulk multifractality spectra 
for p — are compared with the analytical asymptotics, 
Eqs. (|25|) . (|27|) . supplemented by the symmetry relation 
([2"5]) . Again, a very good agreement is seen, except for a 
vicinity of q = 1/2, where Eqs. ([25]), ([27]) break down^ 



V. CONCLUSIONS 

In summary, we have studied the boundary multifrac- 
tality of wave functions in the PRBM model describing a 
critical ID system with long-range hopping. Our findings 
strongly corroborate the ubiquity of the notion of bound- 
ary mutlifractality (recently introduced in Ref. fla ) in the 
context of disordered electronic systems at criticality. We 
have demonstrated, both analytically and numerically, 
that the surface multifractal exponents r|(6,p) are not 
only different from their bulk counterparts, t^(6), but 
also depend on an additional parameter p characterizing 
the reflection of the particle at the boundary. This pe- 
culiarity of the PRBM model is intimately related to the 
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existence of the line of fixed points (labelled by b) in the 
bulk model. Indeed, the freedom in the choice of the am- 
plitude of the boundary "hopping with reflection" term 
is of the same origin as the freedom in the amplitude of 
the power-law hopping in the bulk. 

We close on a somewhat speculative note. The exis- 
tence of a truly marginal coupling (implying a line of 
fixed points) is not unique for the PRBM model. In 
particular, the 2D Dirac fermions in a random vector 
potential^ share this feature. Furthermore, it was con- 
jectured recentl y 28 ' 29 that the quantum Hall transition 
might be described by a particular point on a line of fixed 
points in a related model. Based on our results, it is nat- 
ural to ask whether the emergence of an additional pa- 
rameter p governing the boundary criticality is a general 
property of critical theories with a truly marginal cou- 
pling. The work in this direction is currently underway^ . 
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